Two-dimensional melting far from equilibrium in a granular monolayer 
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We report an experimental investigation of the transition from a hexagonally ordered solid phase 
to a disordered liquid in a monolayer of vibrated spheres. The transition occurs as the intensity of the 
vibration amplitude is increased. Measurements of the density of dislocations and the positional and 
orientational correlation functions show evidence for a dislocation-mediated continuous transition 
from a solid phase with long-range order to a liquid with only short-range order. The results show 
a strong similarity to simulations of melting of hard disks in equilibrium, despite the fact that the 
granular monolayer is far from equilibrium due to the effects of interparticle dissipation and the 
vibrational forcing. 
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The nature of the melting transition of a two- 
dimensional solid has been extensively studied since the 
pioneering work of Kosterlitz and Thouless[l| and subse- 
quent extensions by Halperin, Nelson, and Young@,l3, 
known collectively as the KTHNY theory. This contin- 
ues to be an active area of research 0, § , and investiga- 
tions of colloidal particles suggest that the specific form 
of interparticle interactions can play a role in determin- 
ing the nature of the melting transition^. The partic- 
ular system of hard disks has been extensively studied 
by computer simulation, and there is considerable evi- 
dence that the melting transition is well described by the 
KTHNY scenario (see e.g. 0,0 and references therein). 
That scenario describes two-dimensional melting as con- 
sisting of two continuous transitions, the first from a crys- 
talline phase characterized by algebraically decaying po- 
sitional order and long-range orientational order to a hex- 
atic phase characterized by short-range positional order 
and algebraically decaying orientational order. The hex- 
atic phase melts into a liquid with only short-range order. 
Both transitions are mediated by topological defects, the 
first by unbinding of dislocation pairs, the second by un- 
binding of disclination pairs (see @ for a review). 

Single layers of identical macroscopic spheres subjected 
to mechanical vibration undergo transitions that share, 
at least superficially, many of the characteristics of equi- 
librium atomic and molecular melting transitions [1(1 fllL 
Il2l IT^| . The external vibration is necessary to keep the 
macroscopic spheres in motion because of the inevitable 
loss of energy to friction and inelastic collisions. As a 
result of the flow of energy into the system from the ex- 
ternal forcing and out of the system through dissipation, 
detailed balance is violated, and the results of equilib- 
rium statistical mechanics are not generally applicable. 
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We have shown directly that the fluctuating velocities of 
the individual spheres^! and their correlations ^ij dif- 
fer significantly from those found in equilibrium systems. 
Thus it may be expected that the dynamics of the melting 
transition observed for driven, inelastic spheres may be 
quite different from those observed in equilibrium. How- 
ever in at least one instance it has been shown that in cer- 
tain model systems the Gibbs distribution and detailed 
balance can be recovered by coarse-graining over small 
length scales [THj. 

In this paper we report an investigation of the melt- 
ing transition of a monolayer of identical spherical parti- 
cles confined between two vibrating plates. As reported 
previouslv|lli Il3[. a granular monolayer at moderately 
high density and moderately low vibration amplitude 
forms a ordered array, with each particle fluctuating 
about sites arranged on a hexagonal lattice that extends 
across the entire cell. As the vibration amplitude is in- 
creased, the ordered lattice 'melts', and the spheres dis- 
play disordered, fluctuating, liquid-like dynamics. Here 
we show that this transition is continuous and is consis- 
tent with the KTHNY melting scenario. In particular, in 
the 'solid' phase at low vibration amplitude the positional 
correlation function decays slowly, consistent with an al- 
gebraic decay, and the orientational correlation function 
appears to be long-ranged. The hexagonal lattice in this 
regime is nearly defect free. As the vibration amplitude 
is increased, the density of defects increases dramatically, 
but isolated disclinations are not observed. In this regime 
the positional correlation function decays more rapidly, 
while the orientational correlation function decays very 
slowly, consistent with an algebraic decay. This phase has 
the characteristics of the hexatic phase observed in equi- 
librium 2D melting. Finally, at still higher amplitudes, 
isolated disclinations appear, and only short ranged order 
is observed. Increasing the vibration amplitude reduces 
the particle density in the observation region, and the 
density at which the transitions occur is very similar to 
the densities of the transitions observed in simulations 
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of equilibrium hard-disk melting. These results suggest 
that there is a deep similarity between the coarse-grained 
dynamics of the far from equilibrium system of vibrated 
spheres, and the equilibrium dynamics of hard sphere 
systems. 

The experimental apparatus has been described 
previously 0, 0|. Briefly, a 20 cm diameter smooth, 
rigid aluminum plate is mounted horizontally and vi- 
brated sinusoidally in the vertical direction by an elec- 
tromagnetic shaker. The plate is carefully leveled, and 
the acceleration is uniform across the plate to better than 
0.5%. The plate is machined to be flat, but in fact has 
a slight bowl shape, the importance of which will be dis- 
cussed below. Spheres of 316 stainless steel with a di- 
ameter of 1.588 ± 0.0032mm are placed on the plate, 
confined on the sides by a 1.15cm wide aluminum rim 
and from above by a antistatic coated Plexiglas lid 2.54 
mm (1.6 ball diameters) above the bottom plate. The 
A high resolution camera (Pulnix TM1040) placed above 
the center of the plate is used to locate particle positions 
with a precision of less than 0.01 ball diameters (0.17 
pixels) by calculating intensity-weighted centers of the 
bright spots produced by reflections of the illuminating 
light. The area imaged is approximately 1/8 the total 
area of the plate. 

We have previously reported a phase diagram for this 
system as a function of frequency and amplitude for dif- 
ferent particle densities without the lid present [l3L Il6j. 
The vibration amplitude is described by the dimension- 
less acceleration V — Acu 2 /g, where A is the amplitude 
of the plate oscillation, u> is the angular frequency of the 
plate oscillation, and g is the acceleration due to gravity. 
For the results described below, the vibration frequency 
is ui/2tt — 90 Hz. We have previously shown that the 
phase diagram in this regime is not extremely sensitive 
to freauencv|l3j|. The average two-dimensional particle 
density, p — N/N max , where N is the total number of 
particles on the plate and N max is the maximum num- 
ber of spheres in a single hexagonal close packed layer 
{Nmax = 11270 for this setup), is p = 0.787 for the re- 
sults reported below. 

At low accelerations, the spheres are arranged in a 
hexagonal lattice, with a single domain extending across 
the entire cell. Fig. 2Ji shows the instantaneous positions 
of the particles at T = 0.85. The particle positions fluc- 
tuate continuously, but no particle rearrangements are 
observed. In order to identify defects in the hexagonal 
array, Voronoi constructions are used to identify nearest 
neighbors [2lJ. The neighborhoods marked in green and 
red in FigQJi represent a 5-fold and a 7- fold coordinated 
defect (disclinations), which together represent a dislo- 
cation. At this acceleration the density of dislocations 
is negligible, and the topological order is essentially per- 
fect. As the amplitude of the acceleration is increased 
above T = 1, the density of topological defects increases 
dramatically. Fig. ^ shows a snapshot at V = 1.1, 
where dislocation pairs are common and the ordering is 
significantly disrupted. At still higher accelerations, iso- 



lated disclinations are present (Fig. ^), signifying the 
presence of a disordered liquid. We do not observe any 
evidence for phase separation, which would be expected 
in a discontinuous phase transition. This is in marked 
contrast to the behavior observed when the system can 
form a second layer . w here a clear liquid-solid coexis- 
tence is observed [l9l l20| . There is a small (less than 1%) 
variation in density across the image, due to the slight 
bowl-shape of the plate. The overall density in the im- 
aged region decreases as the shaking amplitude increases, 
because more spheres are pushed out to the edges of the 
plate out of the field of view of the camera. The decrease 
is approximately linear in T, from p = 0.838 at T = 0.85 
to p = 0.791 at T = 1.25 

Figure [21 shows the density of dislocations (circles) and 
6-fold coordinated sites (diamonds) as a function of ac- 
celeration. The low density of dislocations for T < 1.0 
may be due to finite size effects, and in any case has little 
effect on the length scales that are accessible in our ex- 
periment. For r just over 1.0, the density of topological 
defects begins to increase approximately exponentially, 
doubling about every 5% increase in T. The increase 
slows dramatically above T = 1.15, which is coincident 
with the appearance of isolated disclinations. 

The effect of the topological defects on the positional 
ordering of the spheres can be seen through the pair cor- 
relation function G{r) =< p{r')p{r' + r) > / < p(r') > 2 , 
where p is the particle density and the average is taken 
over the spatial variable r'. Fig. [3] shows the pair cor- 
relation function for T = 0.85, 1.1, and 1.2, the same 
accelerations shown in Fig. [21 At low accelerations, the 
positional correlations decay very slowly, consistent with 
an algebraic decay. As the acceleration is increased, the 
decay becomes considerably more rapid, and by T = 1.2 
the positional order is clearly short-ranged. While un- 
ambiguous discrimination between algebraic and expo- 
nential decays is not possible due to the the limited size 
of the system, the evolution of the positional correlation 
function is consistent with the KTHNY melting scenario. 

The topological defects that are responsible for the de- 
cay of the positional correlations in the hexactic phase 
preserve orientational order, so the orientational correla- 
tions in the hexactic phase are predicted to decay slowly. 
Fig. 0]shows the orientational order correlation function, 
96(r) =< ^(r')^(r' + r) > /g(r), where 

N 1 

k j j 

where the first sum runs over the N particles in each im- 
age, and the second sum runs over the Nj neighbors of the 
kth particle in the image. 6kj is the angle of the bond be- 
tween the fcth particle and its jth neighbor. At T = 0.85, 
the orientational order is essentially constant, indicating 
true long-range orientational order. At T = 1.1, there is 
a noticeable decrease in the correlation function at large 
separations. As the acceleration is increased further, the 
rate of the decay increases dramatically, and by T = 1.2, 
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FIG. 1: Particle positions and topological defects. (a) 
T = 0.85, (b) T = 1.1, (c) T = 1.2 Neighborhoods of 6-fold 
coordinated particles are colored yellow, 5-fold green, 7-fold 
red, 8-fold white, and 4-fold blue. Isolated disclinations are 
present in (c), e.g. the 7-fold coordinated particle near the 
top in the center or the 5-fold coordinated particle slightly 
left of the center. 



only short-range orientational order is observed, indicat- 
ing an isotropic liquid. 

The correlation functions in the transition region are 
better fit by an algebraic decay than an exponential, 
although an exponential decay with a long correlation 
length cannot be ruled out. The open squares plotted 
in Fig. [21 show the result of a fit gs(r) = g6(0)r^ v . Be- 
low r = 1 the exponent exponent is negligible, indicating 
long range orientational order. The exponent increases as 
the transition to an isotropic liquid is approached, as in 
the KTHNY melting scenario. In the equilibrium theory, 




0.8 0.85 0.9 0.95 1 1.05 1.1 1.15 1.2 1.25 



FIG. 2: Left axis: Density of dislocations (circles) and 6-fold 
coordinated sites (diamonds) as a function of acceleration. 
Right axis: Exponent from a fit of the orientational order 
correlation function to an algebraic decay, g&(r) oc r~ v (open 
squares) . 
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FIG. 3: Pair correlation function for (a) F = 0.85, (b) T = 
1.1, and (c) r = 1.1. 

77 reaches a maximum of 1/4 at the melting transition, 
at which point ge{r) decays exponentially. In our data 
fits producing values of 77 above 1 indicate short range 
order (see Fig. 0Ji), so within the limitations of our sys- 
tem size and uniformity, these results are consistent with 
the KTHNY melting scenario. 

The close correspondence between the melting transi- 
tion we have observed and the behavior of equilibrium 
systems, combined with our recent observation of a first- 
order-like transition to a two-layer solid^J l2(i| . is a 
strong indication that some aspects of the theory of equi- 
librium phase transitions should be applicable to this far- 
from-equilibrium system. Both of these results are closely 
analogous to well understood transitions in equilibrium 
confined colloids, despite the presence of strong non- 
cquilbrium effects in the velocity distributions ll3l Il7| . 
velocity correlations |14). and energy eauipartition|l9| . al- 
though the results on the two-layer solid also show that 
differing dissipation rates in coexisting phases has a sig- 
nificant effect on the phase diagram. The incorporation 
of these effects into a predictive statistical mechanics of 
nonequilibrium phase transitions remains an open chal- 
lenge. 

The fact that the melting is driven by increasing the 
amplitude of the vibration, and therefore the average ki- 
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FIG. 4: Orientational order correlation function for (a) V = 
0.85, (b) T = 1.1, and (c) F = 1.1. 

netic energy of the spheres, might seem to suggest an 
analogy between the granular temperature (usually de- 
fined as proportional to the mean kinetic energy per par- 
ticle) and the thermodynamic temperature. However the 
stainless steel balls behave as hard spheres, to a very good 
approximation. In equilibrium hard sphere systems, tem- 
perature is not a relevant variable, and phase transitions 
are driven exclusively by density changes. In fact we be- 
lieve a similar mechanism is at work in the results shown 
here. As the shaking amplitude is increased, the effec- 
tive density of the layer changes in two ways. At low 
vibration amplitudes, the balls stay near the plate, and 
the horizontal distance between the centers of the spheres 
cannot be less than the sphere diameter. When the shak- 
ing amplitude is increased, the balls begin to explore all 
of the volume available to them. (As described above, 
the balls are confined from above by a lid that sits 1.6 
ball diameters above the bottom plate.) The horizontal 
distance between the centers of spheres that are not in 
the same horizontal plane can thus be less than the ball 
diameter, effectively reducing the 2D area of each ball, or 



equivalently decreasing the normalized density. Another, 
likely more important effect, arises from the slight bowl 
shape of the plate. The balls tend to settle towards the 
center, but the amount of settling is reduced as the accel- 
eration is increased. Our data was acquired in the center 
portion of the plate, and shows a decreasing density as 
the amplitude is increased. The local density at which 
the transition is observed and the width of the hexatic 
region are in fact very close to those observed in simu- 
lations of hard disk systems^- In an ideal system, the 
plate non-uniformity could be eliminated, but the ball- 
overlap effect would persist. 

Careful measurements of two-dimensional melting in 
magnetic bubble arrays have also observed a continuous 
hexatic-to- liquid melting transition^]- Like the gran- 
ular layer described here, the magnetic bubble arrays 
are far from equilibrium, in the sense that true ther- 
mal fluctuations are irrelevant. In that system, however, 
each magnetic bubble is forced 'microscopically', and ex- 
hibits behavior that is indistinguishable from Brownian 
motion. Thus it is expected to be a good model of equi- 
librium behavior. In the granular layer, however, the de- 
tails of the forcing have a strong influence on the micro- 
scopic dynamics |l4j, which typically differ significantly 
from those observed in equilibrium svstems|l3l 1141 Il6 | . 
The results presented here show that despite these de- 
viations from equilibrium systems, the 'macroscopic' be- 
havior of the granular layer follows the equilibrium melt- 
ing scenario, which suggests that there may be a coarse- 
grained functional that could play the role of the equilib- 
rium free energy in a quantitative theory of the observed 
phase transition. 
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